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ABSTRACT 

The source counts of galaxies discovered at sub-millimetre and millimetre wavelengths 
provide important information on the evolution of infrared-bright galaxies. We com- 
bine the data from six blank- field surveys carried out at 1.1 mm with AzTEC, totalling 
1.6 deg 2 in area with root-mean-square depths ranging from 0.4 to 1.7mJy, and de- 
rive the strongest constraints to date on the 1.1mm source counts at flux densities 
^noo = 1 — 12mJy. Using additional data from the AzTEC Cluster Environment 
Survey to extend the counts to Shoo ~ 20mJy, we see tentative evidence for an en- 
hancement relative to the exponential drop in the counts at Shoo ~ 13mJy and a 
smooth connection to the bright source counts at > 20 mJy measured by the South 
Pole Telescope; this excess may be due to strong lensing effects. We compare these 
counts to predictions from several semi- analytical and phenomenological models and 
find that for most the agreement is quite good at flux densities > 4mJy; however, 
we find significant discrepancies (> 3c) between the models and the observed 1.1 mm 
counts at lower flux densities, and none of them are consistent with the observed 
turnover in the Euclidean- normalised counts at Shoo ^ 2mJy. Our new results there- 
fore may require modifications to existing evolutionary models for low luminosity 
galaxies. Alternatively, the discrepancy between the measured counts at the faint end 
and predictions from phenomenological models could arise from limited knowledge of 
the spectral energy distributions of faint galaxies in the local Universe. 

Key words: sub-millimetre: galaxies - galaxies: starburst - galaxies: evolution - 
galaxies: high-redshift 
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1 INTRODUCTION 

Understanding how star formation evolved over the his- 
tory of the Universe is one of the main goals of extra- 
galactic astronomy today. Dust-obscured star formation is 
known to be a major contributor to the cosmic star for- 
mation history, as the cosmic infrared background (CIRB) 
accounts for ~ 50 per cent of the total extragalactic back- 
ground light (|Puget et al.l [1996). Galaxies that are se- 
lected by their redshifted, thermal dust emission at sub- 
millimetre ( submm) and millimetre (mm) wavelengths, here- 
after SMGs JSmail. Ivison fc Blainlll997J ; lHughes et al.ll998l ; 



iBarger et al.lll998T ). are therefore thought to play a major 
role in the rapid build-up of the stellar populations within 
massive systems. 

SMGs are predominantly high-redshift (z > 1), dust- 
obscured galaxies whose far-infrared (FIR) luminosities 
(£fir > 10 12 Lq) imply high star formation rates (SFRs) of 
> 100 Mq yr _1 ; it is therefore generally believed that SMGs 
are observed during an important starburst or active phase 
in their evolution, en route to be coming massive e lliptical 
galaxies at z = (see review by iBlain et alj 120021 ). Since 
the rest-frame peak of the spectral energy distribution at 
A ~ 100 /mi is increasingly redshifted into the submm/mm 
observing bands with increasing distance, there is a strong 
negative fc-correction for surveys carried out at A > 500 /im. 
SMGs over a wide range in redshift (1 < z < 10) are thus 
readily detected in deep, wide area surveys at these wave- 
lengths, and consequently, their basic statistical properties 
- such as their number density, redshift distribution, and 
clustering strength - hold important clues to how the most 
massive galaxies assemble over time. 

The source counts of SMGs as a function of flux den- 
sity provide strong constraints for modelling the formation 
and evolution of IR-bright galaxies. There are two differ- 
ent methods for incorporating such constraints into mod- 
els. The first method, often referred to as semi-analytical 
or forward evolution models, typically uses numerical simu- 
lations to describe the gravitational collapse of dark mat- 
ter, combined with semi-analytical recipes to govern the 
evolution of baryonic processes within a galaxy, with some 
models including complex processes such as feedback from 
supernovae (SN) and/or active galactic nuclei (AGN) (e.; 
iGranato et "ail 12004 iBaugh et alll2005l : lLacev et al1l201C 
The second type are phenomenological models - or para- 
metric backward evolution models - which make use of 
observational constraints (such as source counts and red- 
shift distributions for galaxy populations selected at dif- 
ferent wavelengths) to derive a model for the evolution 
of the luminosity function of galaxies, considering differ- 
ent populations of galaxies and spectra l energy distribu 



cnt populations 01 gala xies and spectra l energy distribu- 
tions (SEDs) (e.g. [ Valiante et al. 20091; Rowan-Ro binsonl 



120091 ; iPearson fc Khanl 120091 ; iBethermin et all l201lT ). The 
source counts of SMGs measured from both long wave- 
length (A = 850 — 2000 /im) ground-based surveys (e.g . 



Coppin et al. 2006; Bertoldi et ail 120071; IWeifi et all [20091: 



Austermann et al.l |2010| ; IVieira et al.l 120101 ; iMarriage et al.l 
20111 ) and shorter wavelength (A ^ 500 /i m) surveys from 



balloon- or space-based observatories (e.g. Patanchon et al.l 
l2009l ; IClements et al.ll2010l ; IOliver et al.l2010r i require strong 
evolution in the properties of IR-bright galaxies. 

In many ways, surveys at 24 — 500 /mi from the Spitzer 



Space Telescope and the Herschel Space Observatory have 
surpassed those at longer wavelengths in terms of statisti- 
cal power. The large areas combined with the significant 
depths of these surveys make them sensitive to a much 
broader range of galaxy types, whereas existing surveys at 
longer wavelengths are limited to only the most luminous 
systems. Howeve r, due to well known selection effects (e.g. 
IBlain et al.l 120021 ). longer- wavelength data probe, on aver- 
age, higher-redshift and/or colder galaxies, both of which 
are important componen ts for understanding galaxy evolu- 
tion. Several groups (e.g. iDevlin et al.ll2009l ; ICharv fc Popel 
120101 ) have shown that Spitzer/MlPS 24 /im selected galax- 
ies, which are predominately at z < 1.5, account for 55- 
95 per cent of the CIRB at 70-500 fim. However, these 
sources ac count for only ~ 30 per cent of th e CIRB at 
A = 1mm (|Scott et al.ll2010l : iPenner et al.ll201lF ). This sug- 
gests that galaxies selected at 24 /jm, even in the deepest 
surveys, largely miss dust-obscured star formation activity 
at z > 1.5. The study of SMGs at A ^ 850 /mi is thus es- 
sential to improving our understanding of the bulk of star 
formation taking place at higher redshifts. 

There have been several deep (in some cases, confusion- 
limited) surveys carried out at 1.1mm with AzTEC on the 
James Clerk Maxwell Telescope (JCMT) and the Atacama 
Submillime ter Telescope Experiment (ASTE; lEzawa et al.l 
2004, 2008). In this paper, we combine all previously pub- 
lished blank-field survey data taken with AzTEC in order to 
determine the strongest constraints to date on the number 
density of SMGs detected at A Js 850 /im. Given the large to- 
tal area of these combined fields, the uncertainty in the mea- 
sured source counts from cosmic variance is very low. The 
1.1 mm source counts presented here thus provide important 
information on the highest-redshift, IR-bright galaxies, and 
can be used for improving models of galaxy evolution. 

This paper is organised as follows: in Section [2l we pro- 
vide a summary of the blank-field surveys used to derive 
the combined-field source counts; in Section [31 we describe 
the bootstrap sampling method used to derive the 1.1mm 
source counts, and we discuss how we incorporate systematic 
uncertainties from cosmic variance and flux calibration into 
our total error estimates; in Section [4j we discuss estimates 
of the source counts at very high flux densities determined 
from the South Pole Telescope (SPT) and other AzTEC sur- 
veys; we compare the combined- field 1.1 mm source counts 
with predictions from current galaxy evolution models in 
Section [5] we discuss these results in Section [5] and we sum- 
marise this work in Section [7] 



2 SUMMARY OF AZTEC BLANK-FIELD 
SURVEYS 

We select the six individual blan k-field surveys carried out 
with AzTEC (|Wilson et al.ll2008l ) on the JCMT and ASTE 
from 2005 to 2008; each is briefly described below. Table [l] 
lists these fields, the telescope used, map area, depth, and 
number of SMGs detected in each. Note that we do not 
use AzTEC surveys of fields towards known over-densities, 
such as the AzTEC/ASTE map of the SSA-22 fie ld towards 
a proto-cluster at z = 3.1 JTamura et al.l 120091 ) . and the 
AzT EC/JCMT map of the MS-0451. 6-0305 cluster at z = 
0.54 (|Wardlow et al.l l2010f ). Our intent is to produce the 
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strongest constraints on the unbiased 1.1 mm source counts, 
extracted from "blank" fields with no prior known over- or 
under-densities. 

AzTEC map sensitivities tend to decrease from the map 
centre to the edges due to the scanning strategies typically 
employed. For uniformity, we consider the "50 per cent cov- 
erage region" for all fields, which encompasses all pixels in a 
map for which the coverage (i.e. the inverse variance weight) 
is ^ 50 per cent of the maximum coverage. This ensures that 
we are only considering regions of the maps that are well 
sampled by several detectors, with good cross-linking, and 
where the noise properties are uniform. The area and the 
range of la root-mean-square (rms) depth for each field are 
listed in Table [T] The combined fields result in a total area 
of 1.6 deg 2 mapped to la — 0.4-1.7 mjy beam -1 . The res- 
olutions of the JCMT and ASTE data at A = 1.1mm are 
6 = 18" and 28" (full width at half maximum), respectively. 

All of the AzTEC data were reduced using the stan- 
dard customised d ata reduction pipe line in IDL; this is de- 
scribed in detail in lScott et al.l (J2008J), so we do not describe 
it here. We recently derived an improvement in our esti- 
mated transfe r func tion for point sources, as discussed in 
iDownes et al.l (J2011T ). This typically results in an increase 
of 10 — 30 per cent in the measured flux densities of point 
sources detected in the maps. This also increases the noise 
in the maps by roughly the same amount, such that the 
number of SMGs detected based on a signal-to-noise (S/N ) 
threshold does not change significantly. IDownes et al.l (J201l|) 
provide revised source lists for the majority of previously 
published catalogues. This will also result in a shift of the 
source counts published prior to this correction, which af- 
fects all of the fields considered here with the exception of 
COSMOS (|Aretxagaet alj|201ll ); however, the effect on the 
source counts is smaller, since the hi gher noise le vels make 
flux boosting effects stronger (see Austcr mann et al.l 2009, 
|2010j), and in turn, the de-boosting corrections are larger. 
AzTEC/ JCMT survey of the GOODS-N field The AzTEC 
survey of the Great Observatories Origins Deep North 
(GOODS-N) field was carried out during the 2005 to 2006 
observing ca mpaign on the JCMT and is presented in 
iPereraet all ((2008). GOODS-N is one of the most studied 
fields at all wavelengths, and much work has been done to 
identify the multi- wavelength co unterparts to the SM Gs dis- 
covered in the AzTEC survey |Chapin et al.l 120091 ). These 
data have therefore been used extensively to characterise 
the reds hift distribution, AGN fraction, etc. of mm-selected 
sources |Yun et al.ll201ll . Johnson et al. 2012). The revised 
catalogue for this survey, which incorporates the improved 
transf er function estimate, is presented in IDownes et al.l 
J2011I ). 

AzTEC/ JCMT survey of the LH field A region in the 
Lockman Hole (LH) field was observed by AzTEC on 
the JCMT during the 2005 to 2006 observing campaign 
as part of the SCUBA HAlf-Degree Ex tragalactic Survey 

S SHADES) project and is described in lAustermann et al.l 
20101 ). SHADES consists of two discontiguous fields: the 
LH, and the Subaru/XMM-Newton Deep Field (SXDF). We 
do not use the AzTEC/ JCMT map of the SXDF in our com- 
bined source counts analysis, since it largely overlaps with 
the ASTE survey of the same field (see below) and is con- 
siderably shallower. As with G OODS-N, the revised source 
catalogue for LH is presented in IDownes et al.l |201lr ). 



AzTEC/ ASTE survey of the GOODS-S field The AzTEC 
survey of the GOODS-South (GOODS-S) field was carried 
out on the ASTE telescope du r ing th e 2007 observing run 
and is presented in IScott et "all (|2010h . As with GOODS-N, 
extensive efforts to identify multi-wavelength counterparts 
for AzTEC/GOODS-S SMGs, and to derive the redshift dis- 
tribution, SFR, and stellar mass propertie s of SMGs from 
these data have already been carried out (|Yun et al.ll201ll . 
Johnson et al. 2012). The revised ca talogue using the new 
transfer function estimate is given in IDownes et al.l (|2011T ) ■ 
In general, the fractional increase in the measured source 
flux densities and map noise is lower for the ASTE maps 
compared to the JCMT maps. 

AzTEC/ ASTE survey of the ADF-S The AzTEC map of 
the Akari Deep Field South (ADF-S) was built up over the 

2007 and 2008 observin g runs on ASTE and is discussed in 
lHatsukade et al.l (|201lh ■ This is the deepest map used in our 
analysis and therefore puts strong constraints on the faint 
end of the source counts. 

AzTEC/ ASTE survey of the SXDF The AzTEC survey 
of the SXDF field carried out on ASTE during 2007 and 

2008 is a slightly smaller but considerably deeper survey 
than the AzTEC SXDF map taken as part of the SHADES 
project on the JCMT. The source cou nts from this survey 
are presented in lHatsukade et al.l (|201lh . and, similar to the 
ADF-S, these data provide strong constraints on the faint 
end of the source counts. 

AzTEC/ ASTE survey of the COSMOS field The largest 
survey used in our combined source counts analysis is the 
AzTEC survey of the COSMOS field, carried out during th e 
2008 observing campaign on ASTE (|Aretxaga et al.llioilT ). 
This survey almost completely encompasses the smaller, 
shallo wer AzTEC map taken with the JCMT in 2005 to 
2006 (jScott et al.ll2008T ). so we do not use the latter in our 
analysis. Being a factor of > 2 larger than the other sur- 
veys considered here, this field provides the strongest con- 
straints on the bright end of the 1.1 mm source counts. The 
m ap used in our a n alysis here is the same as that presented 
in lAretxaga et al.l (I201lr). which us ed the improved trans- 
fer functionoflDwnes^ai] (|201lh . Like the two GOODS 
fields, COSMOS is one of the best studied regions at all 
wavelengths, and the SMGs detected in this field have been 
used to study the properties of the SMG population at com- 
plementary wavelengths (Johnson et al. 2012). 



3 1.1 MM SOURCE COUNTS FROM 
COMBINED BLANK-FIELDS 

3.1 Bootstrap sampling method 

To derive the 1.1mm source counts, we adopt the stan- 
dard bootstrap sampling method that has been used exten- 
sively in the past for extracting the cou nts from AzTEC sur - 
veys. This method, first introduced bv lCoppin et al.l (|2006h 
and further developed for use with AzTEC d ata, is de- 
scribed in great detail in lAustermann et al.l (J2009J) and 
lAustermann et al.l (|2010j), and we briefly summarise it here. 
Using the source catalogue from one or more surveys 
and assuming a prior distribution for the source counts based 
on the best-fit Schech ter function to the COSMOS counts 
l|Aretxaga et alj|201lf ). we construct posterior flux distribu- 
tions (PFDs) for each source that are sampled at random 



4 Scott et al. 

Table 1. Summary of AzTEC blank-field surveys. The columns are: 1) the field name; 2) the telescope used for the survey; 3) the area 
of the survey within the 50 per cent uniform coverage region (see Section [2]l ; 4) the range of rms noise within the 50 per cent uniform 
coverage region; 5) the number of SMGs detected in the 50 per cent uniform coverage region whose probability of de-boosting to < mjy 
is P(S < 0) ^ 0.20 for sources detected in JCMT surveys and P( S < 0) ^ 0.05 f or sou rces detected in ASTE surveys (see Section [3~H ; 6) 
the root cosmic variance for the survey, following the estimate of lMoster et alj J201ll . Section 13.311 ; and 7) the previous paper describing 
the survey. 
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0.4- 1.7 
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in order to determine intrinsic flux densities for the sources 
in the catalogue; these are then binned to derive the differ- 
ential and integrated source counts. Only sources that pass 
the "null threshold" test are sampled in order to avoid in- 
cluding a large number of false positives; that is, we only 
sample sources for which the probability that their intrinsic 
flux is less than zero is P(S < 0) ^ 0.20 for sources de- 
tected in JCMT surveys and P(S < 0) ^ 0.05 for sources 
detected in ASTE surveys. The more stringent limit for the 
ASTE data is imposed due to larger systematic s from con- 
fusio n in estimating the PFDs (see discussion in lScott et al] 
120101 ). This process is repeated 20,000 times in order to suf- 
ficiently sample the source count probability distribution, 
and the mean and 68.3 per cent confidence interval for the 
counts in even-spaced, 1 mjy-wide flux bins are computed 
from these iterations, giving the raw source counts. These 
raw counts are then corrected for incompleteness, which is 
estimated through simulation by calculating the recovery 
rate as a function of flux density for simulated sources in- 
jected (one at a time) into the map (see lScott et al.ll2010l . for 
example). We use the same source detection algorithm and 
null threshold test on these simulated sources as that used 
for the real catalogues to quantify the survey completeness. 
These corrected counts are then divided by the survey area 
to determine the differential (dN/dS) and integral (N(> S)) 
source counts. 

We compute the counts only for flux densities ^ 1 mjy; 
at lower flux densities, completeness is too low (< 10 per 
cent) and difficult to estimate. While the completeness in 
our three shallowest fields (GOODS-N, LH, and COSMOS) 
is < 5 per cent at lmjy, the deeper surveys are 15-30 per 
cent complete at this same flux level, ensuring that the low- 
flux end of the counts will not be subject to significant errors 
from biases in our completeness estimate. 

Another common method for extracting source counts 
from this type of low-resolution, confusion-limited sur- 
vey is the probability of deflection, or u P(Dy , approach, 
and it has been commonly employed for extracting counts 
from recent BLAST and Herschel-SPI KE surveys (e.g. 
IPatanchon et al.ll2009l ; Idenn et aI.ll201Ch . The P(D) tech- 
nique avoids certain biases inherent in the bootstrap sam- 
pling method - namely, the bias to the counts from the as- 
sumed prior distribution, and biases from the assumption 



that each detected "source" really represents the emission 
from a single galaxy. Also, in principle, the P(D) method 
allows an estimate of the source counts at fainter flux den- 
sities, below the detection limit of individual point sources. 
On the other hand, source counts determined from the P(D) 
approach must use piecewise models, where the differential 
counts at selected "nodes" (i.e. fixed flux densities) are free 
parameters, and the nodes are connected by some smooth 
function. Such models may adequately reproduce the ob- 
served fluctuations in a map; however, they are not at all 
physically motivated. While increasing the resolution be- 
tween nodes can reduce the model dependency of the counts, 
this increases the number of free parameters as well as the 
correlations between them. In practice, most groups limit 
the number of nodes so that the fitted parameters are largely 
uncorrelated, at the expense of making their results more 
model dependent; consequently, the formal errors on the fit- 
ted parameters may not always represe nt the true uncer- 
tainty in the cou nts (see discussions in IScott et all l20ld ; 
Idenn et al.ll2010l . for example). Furthermore, while the im- 
plementation of the P{D) method is relatively straightfor- 
ward when the transfer function for point sources is lin- 
ear, this is not the c ase for our PCA-cleaned AzTEC maps 
IjDownes et al.ll201ll ). and the P(D) method thus becomes 
computationally expensive for our data. This is why we have 
chosen to use the bootstrap sampling method instead. 

lAustermann et al] (|2010r i demonstrated that for flux 
bins that are well sampled, the assumed prior used in the 
bootstrap sampling approach is quite weak. Furthermore, 
they show that biases to the counts for sparsely sampled 
flux bins can be effectively removed by an iterative process 
in which the counts extracted from the first pass are used 
to update the prior and PFDs for the source catalogue(s), 
and the bootstrap sampling method is repeated. We use this 
iterative process to decrease the effects of the prior on the ex- 
tracted source counts. The counts from the combined blank 
fields change by less than 2 per cent in all flux bins after only 
three iterations. 

We list in Table [l] the total number of sources in each 
field that were used to extract the counts. The counts from 
the combined blank-fields are shown in Fig. [T] (filled circles) , 
and are listed in Table [2] where the upper and lower error 
bars indicate 68.3 per cent confidence intervals. We give the 
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Table 2. 1.1mm source counts derived from the combined six 
AzTEC blank-field surveys. The first two columns show the flux 
bin centres and corresponding differential source counts, while the 
last two columns show the flux bin minima and cumulative counts. 
The first set of upper and lower errors shown on the counts indi- 
cate the 68.3 per cent confidence intervals considering only statis- 
tical errors (Section l3.1l l. The second set of errors in parentheses 
shows the 68.3 per cent confidence intervals when including sys- 
tematic uncertainties from cosmic variance (Section 13.31 1 and flux 
calibration (Section 13. 4jl , The bright source counts derived from 
the AzTEC Cluster Environment Survey (ACES) are also listed 
(see Section B) . 
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bles[3]and[4] respectively, and we list the standard deviation 
on the counts, (TdN/ds an d cjv(>s) i f° r each bin as well, from 
which the covariance matrices can be determined. We have 
made the source counts and correlation matrices available 
online for public use(j Since the counts are determined from 
bootstrapping off the PFDs of the sources, adjacent flux bins 
are strongly correlated. It is therefore important to use the 
covariance matrix in model fitting; e.g., for fitting a model 
prediction m to observed counts d with covariance matrix 
C, the \ 2 metric is given by \ 2 — (d — m)C _1 (d — m) T . 



3.2 Effects of confusion on the extracted source 
counts 

Usin g the standard definiti on of one source per 30 beams 
(e.g. iTakeuchi fc Ishiil 120041 ). the confusion limit for these 
surveys is Sum = 1.4 and 2.4 mjy for those carried out on 
the JCMT and ASTE, respectively. These correspond to the 
two lowest flux bins in our source counts estimate. In this 
section, we explore potential biases to the extracted counts 
due to confusion effects. 

1 http://www.astro.umass.edu/AzTEC/Scott2012_nc/ 
aztcc_combined_counts_2012.html 



For this purpose, we make fully simulated data-sets for 
each of the six blank fields. These simulated maps have the 
same noise properties as the real data, and are all populated 
with the same source distribution described by a Schechter 
function: 



dN 
dS 



= N 3 



mjy 



3mJy 



-(S-3mJy)/S' 



(1) 



We use values of -V_ m j y = 230 mjy -1 deg -2 , S" = 1.7 mjy, 
and a — —2, which provide a good fit to the observed blank- 
field source counts derived in this paper. Since we are only 
looking to examine general trends in potential biases to the 
counts from source confusion, a priori knowledge of the true 
distribution of the underlying counts is not necessary. We 
populate each simulated map with sources down to a flux 
density limit of 0.1 mjy, where the cumulative counts reach 
> 1 source per beam for both JCMT and ASTE surveys. 
Previous studies (e.g. IScott et al.ll2010l ) have demonstrated 
that the choice of this lower flux density cutoff is not too 
critical, so long as it corresponds to where the integrated 
source counts are > 1 source per beam, since adding fainter 
sources at that point would not change the flux distribution 
in the map. Simulated sources are placed at random posi- 
tions drawn from a uniform spatial distribution. We make 
100 simulated maps for each field and use the bootstrap sam- 
pling method described in Section 13.11 to derive the source 
counts for each of them. 

The results from these simulations are presented in 
Fig. [2] which shows the Euclidean-normalised differential 
source counts averaged over the 100 simulated maps for each 
field separately (top panel) . The counts have been scaled ar- 
bitrarily for clarity, with dotted curves indicating the input 
source distribution from Equation [T] There is evidence from 
these simulations that confusion introduces biases to the ob- 
served source counts. This is more evident in the bottom 
panel of Fig.[5J which shows the fractional difference between 
the input and output source counts. In general, the counts at 
5*1100 > Slim (the confusion limit) are overestimated by ~ 5- 
30 per cent, while the counts at Shoo < Slim are underesti- 
mated by ~ 10-20 per cent. However, these biases are small 
compared to the statistical errors on the derived counts for 
these individual fields; the error bars in the top panel of 
Fig. [2] represent the typical 68.3 per cent confidence intervals 
on the extracted counts for a single simulated map. Con- 
sidering all simulated data-sets, the extracted counts agree 
with the input counts within their 2<r errors > 85 per cent 
of the time, with the exception of the COSMOS simulated 
fields, where the extracted counts at Shoo = 1.4 mjy are 
always significantly underestimated. This large discrepancy 
between the input and output source counts in the lowest 
flux density bin for COSMOS is most likely due to the low - 
and therefore poorly measured - completeness at that flux 
density. 

We next use the simulated maps for each field to make 
100 realisations of the extracted source counts from the six 
fields combined. For each realisation, we randomly select six 
simulated maps, one from each field, and carry out the joint 
bootstrap sampling extraction. These results are also shown 
in Fig. [2] and as expected, we see similar biases to the output 
source counts as seen in each individual field. This bias is 
smallest (4 per cent) for the 2.4 mjy flux density bin, which 
corresponds to the confusion limit for the ASTE surveys. For 
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Table 3. Correlation matrix for the differential counts derived from the combined AzTEC blank-field surveys. The first set shows 
the bin-to-bin correlations when considering only statistical errors from the bootstrap sampling method (Section 13. l|l ; these represent 
the actual correlations between bins in the data themselves. The second, third, and fourth sets show the correlations when systematic 

uncertainties from cosmic variance (Section 13.31 1. flux calibration (Section 13.41 1. and both are included, respectively. The last column in 
all four cases shows the standard deviation on the differential counts for each flux bin, and can be used to compute the covariance matrix 
from these correlations. 

Flux Density VdN/dS 

(mjy) 1.4 2.4 3.4 4.4 5.4 6.4 7.4 8.4 9.4 10.4 11.4 (mjy- 1 deg" 2 ) 

Statistical errors only 

1.4 1.00 66 

2.4 0.61 1.00 27 

3.4 0.21 0.62 1.00 14 

4.4 0.11 0.29 0.67 1.00 8.2 

5.4 0.06 0.15 0.34 0.71 1.00 5.0 

6.4 0.02 0.06 0.15 0.33 0.67 1.00 3.2 

7.4 0.00 0.01 0.05 0.13 0.32 0.72 1.00 2.1 

8.4 -0.01 -0.01 0.00 0.04 0.12 0.32 0.76 1.00 1.6 

9.4 0.00 -0.01 -0.01 0.01 0.05 0.17 0.49 0.81 1.00 1.2 

10.4 0.00 -0.01 -0.02 -0.01 0.01 0.06 0.20 0.42 0.83 1.00 0.87 

11.4 0.00 -0.01 -0.02 -0.01 -0.01 0.01 0.07 0.20 0.54 0.85 1.00 0.70 

Including systematic uncertainties from cosmic variance 

1.4 1.00 80 

2.4 0.80 1.00 31 

3.4 0.57 0.79 1.00 15 

4.4 0.48 0.58 0.79 1.00 8.8 

5.4 0.38 0.44 0.54 0.79 1.00 5.2 

6.4 0.27 0.30 0.35 0.47 0.72 1.00 3.2 

7.4 0.19 0.21 0.22 0.27 0.40 0.73 1.00 2.1 

8.4 0.13 0.14 0.14 0.14 0.19 0.36 0.77 1.00 1.6 

9.4 0.11 0.12 0.11 0.11 0.12 0.22 0.51 0.81 1.00 1.2 

10.4 0.09 0.09 0.08 0.07 0.07 0.09 0.22 0.43 0.83 1.00 0.87 

11.4 0.07 0.08 0.07 0.06 0.05 0.05 0.09 0.22 0.55 0.85 1.00 0.71 

Including systematic uncertainties from flux calibration 

1.4 1.00 94 

2.4 0.58 1.00 27 

3.4 0.20 0.62 1.00 14 

4.4 0.11 0.27 0.67 1.00 9.3 

5.4 0.05 0.14 0.33 0.71 1.00 6.5 

6.4 0.02 0.06 0.16 0.32 0.68 1.00 4.0 

7.4 0.00 0.02 0.06 0.13 0.30 0.67 1.00 2.5 

8.4 0.00 0.00 0.01 0.03 0.11 0.34 0.83 1.00 1.7 

9.4 -0.01 0.00 0.00 0.01 0.04 0.16 0.46 0.75 1.00 1.3 

10.4 0.00 0.01 0.00 0.00 0.00 0.04 0.17 0.42 0.86 1.00 0.92 

11.4 0.01 0.01 0.01 0.00 0.00 0.00 0.05 0.17 0.48 0.79 1.00 0.72 

Including systematic uncertainties from cosmic variance and flux calibration 

1.4 1.00 100 

2.4 0.78 1.00 32 

3.4 0.57 0.80 1.00 16 

4.4 0.47 0.58 0.79 1.00 9.9 

5.4 0.38 0.44 0.55 0.79 1.00 6.7 

6.4 0.26 0.30 0.34 0.45 0.73 1.00 4.0 

7.4 0.18 0.19 0.20 0.24 0.38 0.69 1.00 2.5 

8.4 0.13 0.14 0.14 0.15 0.21 0.40 0.85 1.00 1.7 

9.4 0.12 0.12 0.12 0.12 0.14 0.22 0.46 0.75 1.00 1.3 

10.4 0.09 0.10 0.10 0.09 0.09 0.11 0.20 0.46 0.86 1.00 0.91 

11.4 0.08 0.09 0.08 0.07 0.07 0.06 0.07 0.18 0.44 0.76 1.00 0.72 
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Table 4. Correlation matrix for the cumulative counts derived from the combined AzTEC blank-field surveys. The first set shows 
the bin-to-bin correlations when considering only statistical errors from the bootstrap sampling method (Section 13. l|l ; these represent 
the actual correlations between bins in the data themselves. The second, third, and fourth sets show the correlations when systematic 
uncertainties from cosmic variance (Section 13.31 1, flux calibration (Section 13.41 1. and both are included, respectively. The last column in 
all four cases shows the standard deviation on the cumulative counts for each flux bin, and can be used to compute the covariance matrix 
from these correlations. 
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Figure 1. Differential source counts derived from the six blank-field surveys carried out with AzTEC on the JCMT and ASTE. The 
counts determined from each individual field are as follows: triangles - GOODS-N; squares — LH; diamonds - GOODS-S; inverted 
triangles - ADF-S; circles - SXDF; stars - COSMOS. The counts derived from combining these six fields are shown as black filled circles. 
The counts for the individual blank-fields have been computed in slightly different flux bins for clarity in plotting. All error bars show the 
68.3 per cent confidence intervals determined from the bootstrap sampling method, including uncertainties arising from cosmic variance 
(see Section 13.311 . The uncertainty from flux calibration (Section 13. 4H is not included since all fields were calibrated the same way. The 
shaded region highlights the 68.3 per cent confidence range on the combined counts. The horizontal lines at the bottom right-hand side 
indicate the survey limits for each individual field (dotted) and for the combined counts (dashed). The survey limit corresponds to the 
expected value for which the source counts will Poisson deviate to zero 31.7per cent of the time, given the area of the survey(s). 



the lowest flux bin at 1.4 mjy, the counts are underestimated 
by 17 per cent, and for the bins at Shoo > 3mJy, the counts 
are overestimated by 9-34 per cent. In comparison to the 
statistical errors, the combined-field counts derived from the 
simulated maps agree with the input counts within their 2<r 
errors > 80 per cent of the time. 



the GOODS-N counts are consistently higher than the av- 
erage, while the GOODS-S counts are consistently lower. In 
order to assess the agreement among the counts derived from 
individual AzTEC fields, we must include this uncertainty 
from cosmic variance into the error budget. 



3.3 Cosmic variance 

The source counts determined from each individual survey 
are shown in Fig. [T] Given the limited area surveyed for each 
field, we expect to see variations in the counts from field-to- 
field owing to variations in the underlying large-scale struc- 
ture (also known as "cosmic variance"). Furthermore, given 
the strong bin-to-bin correlations, the counts across all flux 
bins for a given survey vary in the same sense; for example, 



We estimate the expected cosmic variance for each in- 
dividual survey and for the combined blank- fields using the 
prescription described in lMoster et al.l (|201ll ). This estimate 
uses predictions of the underlying structure of cold dark 
matter and the expected bias for a galaxy population - in 
this case, SMGs. This simple recipe depends only on the an- 
gular dimensions of the field (ai, 012), the mean redshift (z), 
redshift bin size (A2), and stellar mass (M*) of the galaxy 
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Figure 2. Results of simulations to test the effects of confusion 
on our source counts extraction, as described in Section l3.2l The 
top panel shows the averaged Euclidean-normalised differential 
source counts extracted from simulated maps for each field and 
for the combined fields, as labeled. For clarity, the counts have 
been offset by a factor of 1000, 100, 10, 0.1, 0.01, and 0.001 for 
the GOODS-N, LH, GOODS-S, ADF-S, SXDF, and COSMOS 
data-sets, respectively. The dotted curves show the model for the 
input source counts (Equation [TJ for comparison. The errors rep- 
resent typical 68.3 per cent confidence intervals on the counts for 
a single simulated map, demonstrating that the biases to the ex- 
tracted counts arising from confusion effects are small compared 
to the statistical errors. The bottom panel shows the fractional 
difference between the input and output source counts, using the 
same symbols as in the top panel. 



population in question. iMoster et al.l (|20Tlf ) have provided 

their software tools for calculating cosmic variance onlinel 2 ] 

This estimate assumes rectangular geometry for the 

survey, which is rarely the case for our fields; however, 

2 http://www.mpia.de/homes/moster/rcsearch 



IMoster et al.l (|201lh show that the geometry makes little 
difference except where the ratio between the short and 
long axes of the survey is < 0.2, which is not the case 
for any of our fields. We therefore assume angular dimen- 
sions for each field, Qi and 0:2, such that the product is 
equal to the area in the 50 per cent uniform coverage re- 
gion, and the ratio approximately matches the geometry of 
the AzTEC map. The mea n redshift, redsh ift bin size, and 
stellar mass are taken from lYun et al.l (|2011f ). who use spec- 
troscopic and photometric redshifts for the SMGs detected 
in both AzTEC/GOODS fields to determine their redshift 
distribution, and stellar mass estimates fro m modelling their 
observed rest-frame UV and optical SEDs. lYun et al.l (|201lh 
find z = 2.5, Az = 1.5, and M* > 10 10 M , where Az and 
the limit on M+ encompass 75 per cent of the SMGs in that 
sample. The root cosmic variance, cr gg , which represents the 
expected fractional error on the counts due to cosmic vari- 
ance, is listed for each field in Table [T] and ranges from 6.5- 
12.9 per cent for the individual AzTEC fields. By combining 
all six fields, totalling 1.6 deg 2 , uncertainties due to cosmic 
variance are reduced to only 3.9 per cent, which is smaller 
than the statistical errors on the counts (^ 6 per cent). 

Since the uncertainty from cosmic variance is com- 
pletely correlated among all flux bins, we cannot simply add 
cr g g in quadrature with the statistical errors on the counts; 
instead, we must consider how including cosmic variance af- 
fects the entire covariance matrix. It is straightforward to 
include this effect within the framework of the bootstrap 
sampling method. For each of the 20,000 iterations, we gen- 
erate a random number drawn from a Gaussian distribution 
with a mean of zero and a standard deviation of <T gg , and 
we apply this fractional deviation to the differential source 
counts uniformly to all bins. The mean, 68.3 per cent con- 
fidence intervals, and the covariance matrix for the counts 
are then computed from the 20,000 iterations in the same 
manner as in the standard bootstrap method described in 
Section 13.11 This way, we broaden the distribution in the 
bootstrapped counts according to the expected degree of 
cosmic variance and can properly trace the effects on the 
bin-to-bin correlations. 

The 68.3 per cent confidence intervals shown by the er- 
ror bars on the differential counts in Fig.[T]include the uncer- 
tainties expected from cosmic variance, for each individual 
field as well as the combined fields. We note that the stan- 
dard deviation of the counts, <J d N/dS (equal to the root of 
the diagonal elements of the covariance matrix) , increases as 



expected {o 2 dN/dS 



T dN/dS 



+ o-g g (dN/dS) ), and the bin- 



to-bin correlations on the counts increase substantially, as 
shown in Tables [3] and U Comparing the counts observed 
in each individual field when both the statistical errors and 
the uncertainties from cosmic variance are included, we find 
that they all agree quite well. 

The m ean redshift and inter-quartile range from 
lYun et al.l (|2011T ) which we use to estimate the cos- 
mic variance agree very well with tho s e der ived from 
other spec troscopic ( Chapm an et al. 20031. 20 05|) and pho- 
tometric (lPopeetal.1 120051 ; lAretxaga et al.l |2003| , 120071 ; 
IWardlow et al.l 120111 ) redshift estimates, which find z — 
2.2 — 2.5 and Az = 1.1 — 1.8. The largest uncertainty in this 
estimate comes from the assumed stellar mass distribution 
of SMGs. Stellar mass estimates for SMGs are highly uncer- 
tain, as they depend strongly on the choice of the stellar syn- 
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thesis model, star-formation history, and initial mass func- 
tion (IMF ) - all of which are not well understood (see dis- 
cussion in lMichalowski et al.ll201.il . for example). The stellar 
mass distribution of AzTEC/GOODS SMGs in lYun et all 
|201ll . see their figure 8) peaks at M* ~ 10 113 M and is 
broadly consistent with the mean stellar masses of SMGs 
estimated in other works, which range from ~ 10 108 to 
10 11 - 8 Mn (iDve et all 120081; IMichalowski, Hiorth fc Watsonl 
|2010| ; iHainline et al.ll201ll ; IWardlow et al.ll201ll) . HoweverT 
the distribution in M* from all of these works is found 
to be quite broad, especia l ly com pared to the stellar mass 
bins used in lMoster et al.l 1)20111 ) for computing the galaxy 
bias; this is why we opt to use a lower stellar m ass limit of 
M*_> 10 10 M Q , which includes 75 per cent of the lYun et all 
IJ20111 ) sample. Increasing the lower limit on the stellar mass 
would increase the expected cosmic variance for these sur- 
veys, as more massive galaxies are more strongly clustered. If 
we instead assume M+ > 10 10,5 M^ (as moti vated to some 
extent by results in IMichalowski et al.l l201ll . see their fig- 
ure 3), we would derive a root cosmic variance of 8.8-17.5 per 
cent for the individual AzTEC surveys and 5.3 per cent for 
the combined fields. 



3.4 Systematic uncertainty from flux calibration 

We must also consider a systematic uncertainty on the de- 
rived source counts arising from uncertainty in the absolute 
flux calibration of our data. For AzTEC data, we determine 
flux conversion factors to convert the raw detector signals to 
flux density units based on several calibration observations 
of p lanets taken oyer a w ide range of atmospheric conditions 
(see lWilson et al.ll2008l . for details). While the random cal- 
ibration error of a n individual observation is ~ 10 per cent 
|Scott et al.l l2010h . co-added AzTEC maps are each built 
from 46 to 325 observations, and the error on our measured 
source flux densities integrates down to 0.3 per cent when all 
observations, and all fields, are considered. However, for all 
of these data we use the same flux calibrators, Uranus and 
Neptune, which have an abs olute uncer tainty on their flux 
densities of a ca i = 5 per cent (|Griffin fc Ortonll 19931 ) . This is 
a systematic uncertainty in the flux scale of our maps that is 
completely correlated among all AzTEC data and therefore 
propagates into the source counts. 

As with cosmic variance, we incorporate this calibra- 
tion uncertainty into the bootstrap sampling method. For 
each of the 20,000 iterations, we generate a random num- 
ber drawn from a Gaussian distribution with a mean of zero 
and a standard deviation of a ca i. We then modify the PFD 
of every source assuming that the observed flux and noise 
change by this fractional amount, consistent with a system- 
atic change in our flux calibration. We use the same method 
of sampling the PFDs as described in Section 13.11 creating 
20,000 realisations of the counts, from which we compute the 
mean, 68.3 per cent confidence intervals, and the covariance 
matrix, as before. 

The correlation matrices for the counts, when includ- 
ing this systematic calibration uncertainty, are shown in 
Tables [3] and [4] The standard deviation on the differential 
counts increases by sj 5 per cent in all flux bins. We find 
that the bin-to-bin correlations on the differential counts do 
not increase substantially; this is because the perturbations 
to the PFDs of the sources that account for the absolute 



calibration uncertainty are small compared to the intrin- 
sic width of the PFDs owing to the low signal-to-noise of 
the detections. It is this latter feature that gives rise to the 
strong correlations seen among the bins before including any 
systematic uncertainties. 

The differential source counts and 68.3 per cent con- 
fidence intervals, when systematic uncertainties from both 
cosmic variance and flux calibration are included, are shown 
in Fig. [3] Tables[21|3] and|3]list the differential and integrated 
source counts, 68.3 per cent confidence intervals, and corre- 
lation matrices, with and without systematic uncertainties. 



4 BRIGHT SOURCE COUNTS FROM ACES 
AND THE SPT 

Over the flux densities for which these AzTEC blank-field 
surveys are sensitive, the source counts are well described 
by a Schechter function (Equation [l]), declining exponen- 
tially with increasing flux density, since highly luminous 
galaxies are quite rare. Recent surveys at submm and mm 
wavelengths covering > 100 deg 2 have been achieved from 
ground- and space-based observator ies, including the At - 
acama Cosmology Telescope (ACT; iMarriage et al.ll201ll ). 
the SPT dVieira et aT1|2010l), and the Hersch el Space Ob- 
servatory (|Eales et al.ll2010l : lOliver et alj|2010r i. These sur- 
veys are turning out large numbers of extremely bright, 
rare objects that are not associated with known nearby 
systems or strong radio sources, having SEDs consistent 
with high-redshift, dusty star-forming galaxies. It has been 
shown that a significant number of these extremely bright 
systems detected by Herschel, with 500 /im flux densities 
5500 > 100 m jy, are strongly lensed by foreground g alaxies 
or structure (|Negrello et al.jhoiol ; IConlev et al.ll201ll ). with 
magnification factors of /1 ~ 10. These lensed galaxies are 
believed to contribute significantly to the source counts at 
flux densities greater than th ose probed by our compara- 
tively small AzTEC surveys. IVieira et al.l (|2010l ) have es- 
timated the counts at S1400 > lOmJy for 1.4 mm sources 
detected by the SPT, excluding those that are associated 
with nearby galaxies discovered by the Infrared Astronomi- 
cal Satellite (IRAS) and those with synchrotron-dominated 
(as opposed to dust-dominated) SEDs. These are shown 
alongside the AzTEC blank-field counts in Fig. [3] where we 
have scaled the SPT counts to 1.1mm assuming a spectral 
index of a = 2.65, which corresponds to the average spectral 
index between the observed flux densities at 1.1 and 1 .4 mm 
for a starburst galaxy at z = 3 (|Yun fc Ca rilli 2002). The 
SPT counts, which are believed to be dominated by strongly 
lensed galaxies, diverge from the exponential fall-off that 
would be extrapolated from the AzTEC blank-field counts. 

During the 2007-2008 observing seasons on ASTE, we 
observed 37 individual, relatively small (~ 300 arcmin 2 ) 
fields centred on known over-dense regions as part of the 
AzTEC Cluster Environment Survey (ACES). This survey 
was designed to study the SMG population towards biased 
fields and includes regions surrounding clusters and proto- 
cluster candidates from z — 0.05 - 6, covering a total area 
of 3.1 deg 2 . A full description of ACES and first results will 
be discussed in Zeballos et al. (2012). 

The ACES fields were not included in our source counts 
estimate since we specifically wanted to avoid known biased 
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regions. Still, each of the ACES maps includes a relatively 
large area, far from the central over-dense region that is 
not expected to be influenced by the cluster environment. 
We see several bright (^ 8mJy) SMGs located far from the 
biased regions of these ACES maps; from these, we estimate 
the bright counts at 1.1mm, which are poorly constrained 
by our combined blank-field surveys. Since our goal is to fill 
in the information on the unbiased source counts in the flux 
range that is not well sampled by the AzTEC blank-field 
counts or the SPT counts, we apply masks and various 
selection criterion to the ACES data. We first exclude any 
SMG detected within 2 arcmin of the centre of the cluster 
or proto-cluster core. At z > 0.3 (the minimum redshift 
of the ACES clusters where 5*1100 ^ 8mJy sources are 
found), 2 arcmin corresponds to ^ 0.5 Mpc, which is larger 
than the expected core radii for massive clusters with to- 



tal masses 10 



M, 







5 COMPARISON WITH GALAXY 
EVOLUTION MODELS 

With these combined AzTEC surveys, we have put the 
strongest constraints to date on the blank-field 1.1mm 
source counts from S1100 = 1 to 12 mjy. These counts pro- 
vide important information for modelling the formation and 
evolution of galaxies. A detailed analysis of how these re- 
sults fit into our current understanding of galaxy evolu- 
tion is beyond the scope of this paper. However, we can 
compare our observed 1.1 mm counts to existing predictions 
from evolutionary models from the literature, many of which 
have used constraints from 1.1 mm source counts in the past 
from smaller and/or shallower surveys, in order to moti- 
vate modifications to existing models in light of these new 
constraints. We provide only a qualitative comparison here, 

test or 



(e.g. iNavarro. Frenk fc White! since common statistical tests (e.g. the Pearson \ 



19951: iLlovd-Davies Ponman fc CannonT l2000l : 

Kav Liddle fc Thomas! l200ll ). This mask will there- 
fore exclude the major sources of potential biases to 
the counts, including those from strong lensing and the 
Sunyaev-Zeldovich effect, as well as mo st cluster members . 
Next, we follow a similar analysis as in I Vieira et al.l (|201y) 
to exclude nearby galaxies and synchrotron-dominated 
sources from this sample. We exclude sources which 
have 2MASS A"-band and/or IRAS counterparts located 
within the AzTEC beam to eliminate low-redshift galaxies. 
We then check for sources associated with bright radio 
objects using the Aust ralia Telescope Co mpact Array 
(ATCA) 20-GHz survey (JMurphv et alJl201Cl ). the Sydney 
University Molonglo Sky S urvey (SUMSS) at 843 MHz 
(|Bock. Large fc Sadler! 1 19991 ). and a unified catalogue of 
radio sources, which combines in formation from FIRST , 
NVSS, WENSS, GB6, and SDSS JKimball fc Ivezidl200Sh . 
None of the Js 8 mjy ACES sources have radio associations, 
implying that their SEDs are likely dominated by dust 
emission. After culling our sample using the above selection 
criterion, we have a total of 35 Shoo ^ 8mJy sources from 
these ACES fields. 

The ACES bright source counts are shown in Fig. [3] 
and are listed in Table [2] Given small sample statistics, 
the errors on the bright source counts derived from the 
ACES fields are large; however, we see tentative evidence 
for a divergence of the counts from an exponential fall-off 



at Six 



> 



13mJy, and a smooth connection between the 



AzTEC blank-field counts at Snoo < 12mJy and the SPT 
counts at Shoo > 20mJy. This upturn may be highlight- 
ing the regime at which the source counts become domi- 
nated by gravitational lensing effects. On the other hand, 
it is possible that the 2 arcmin mask is not large enough to 
remove all cluster members, and the upturn in the counts at 
Shoo > 13mJy is merely reflecting an over-density of SMGs 
in the outer regions of the clusters. The source counts at 
Shoo = 1 — lOmJy derived from the ACES fields, how- 
ever, do not support this; Zeballos et al. (2012) find that 
the ACES source counts are completely consistent with our 
results from blank-fields over these flux densities when the 
inner regions (1.5 arcmin radii) are masked. Since a real over- 
density would be flux-independent, the lack of an excess in 
the ACES counts at Shoo < 10 mjy suggests that the up- 
turn at Shoo > 13 mjy is not due to cluster-member con- 
tamination. 



the Kolmogorov-Smirnov test) do not apply given the strong 
correlations between these binned data. 

Ten such models are shown in Fig. [3] Th ese include pre- 
diction s from the semi-an alytical models of iGranato et al.l 
l|2004h . lBaugh et al.l (|2005h . and lWilman et al.l IJ2010T). and 
from th e phenomenologi c al mo de ls of Rowan- Robinson! 



(2009 ). IPearson fc Khanl d2009h. IValiante et al 



2009), 



2011) 



Franccschi ni et al. (|2010[) . IBethermin et al.l 
Mars den et al.l 1 20111 ). and iRahmati fc van der WerJ 
1 20111 ). All of these models assume the standard ACDM 
cosmology, but with slightly different parameters ranging 
from A = 0.7 to 0.734, Q m = 0.266 to 0.3, and H = 70 to 
75kms _1 Mpc -1 ; we have not scaled the model predictions 
to a uniform cosmology since the differences should be 
minor. The models show a great deal of dispersion among 
themselves, differing by as much as a facto r of four at 
any g iven flux density. With the exception of iBaugh et al] 
(J2005h and IBethermin et al] i|201ll ). all of the models appear 
to be largely consistent (within ~ 3<r) with the observed 
source counts at Shoo > 4mJy, but many are significantly 
discrepant at lower flux densities. 

With the exception of IBethermin et al] l|201ll ) , these 
models did not include the effects of strong lensing on 
the counts (which is expected to become important around 
Shoo > 15mJy), and did not use the bright counts mea- 
sured by the SPT as constraints. The general agreement be- 
tween these models and the ACES and SPT counts is there- 
fore somewhat fortuitous. However, this should not be over- 
interpreted, since the errors on the ACES and SPT counts 
are quite large, and the scaling used to convert the 1.4 mm 
counts from the SPT survey to 1.1 mm is uncertain. We thus 
limit our discussion in this section to the comparison of these 
models to the AzTEC blank-field counts. 

Considering first the semi-analytical models, the model 
that best fits the o bserved 1.1mm source counts at all flux 
densities is that of [Granato et al.l (|2004J), which considers 
the evolution of gas within dark matter halos as driven 
by gravity, radiative cooli ng, and feedback fr om supernovae 
and AGN. Alternatively, IBaugh et all (|2005l ) modelled the 
850 /im counts available at the time within the context of 
hierarchical assembly, where bursts of star formation are 
triggered solely by galaxy mergers; their model requires a 
top-heavy IMF in order to explain the number density of 
SMGs. Although this model over-predicts the 1.1 mm counts 
by > 3cr at Shoo < 5 mjy, minor changes to the model could 
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bring it more in line with our measurements (e.g., changing 
the dust emissivity). 

The semi-analytical model of IWilman et al.l (|2010l ) is 
an extension of their previous work to simulate the ex- 
tragalactic radio continuum sky, including AGN and star- 
forming galaxies, within the fram e work of the ir large-scale 
cluste ring. Unlike iGranato et al.l 1)20041 ) and iBaugh et al.l 
(2005), who use self-consistent radiative transfer calcula- 
tions to describe the absorpti on and re-emission of starlight 
by dust. lWilman et al.l (|2010n use families of SED templates 
and FIR-radio relationships for star-forming and AGN-host 
galaxies to predict the IR-to-submm emission from their 
simulated radio galaxies, similar to what is done in phe- 
nomenological modelling. Their best-fit model accurately 
predicts the Spitzer 24-160 /im counts, as well as the 850 )im 
counts from SCUBA . However, we find that the model of 
IWilman et al.l (J2010f) over-predicts the 1.1mm counts by 
> 3a at 5*1100 < 3mJy. All three of the semi-analytical mod- 
els considered here predate the large Herschel surveys, and 
it would be interesting to see how they fare when compared 
to the counts at shorter submm wavelengths. 

We now compare the observed counts to predictions 
from phenomenological models. The only two that match the 
blank-field counts withi n 3a at all flux densities are those of 
I Rowan- Robinson! |2009|, wi th a redshift beyond wh ich evolu- 
tion is zero of Zf = 4) and P earson fc Khanl (|2009r ) ■ Both of 
these models consider only four populations of galaxies, each 
represented by a single SED, and are therefore some of the 
simplest models presented here. However, they both used the 
publish ed 1.1 mm source c ounts from the AzTEC/GOODS- 
N field |Perera et al.ll2008l ) to constrain their models, so the 
agreement with the combined-field source counts is not sur- 
prising. 

Like those 
iPearson fc Khan 
of IValiante et al 



of iRowan- Robinson! (|2009l ) and 
(2009), the phenomenolo gi cal m odels 
(l2009h . IFranceschini et al.l CJ2010I ), and 



iMarsden et al.l 1 20111 ) used constraints from source counts 



ranging from 24 to 1100 /im, but did not use constraints 
on the 250-500 /xm coun ts fr om Herschel surveys - though 
IFranceschini et all (|2010l ) and lMarsden et all (1201 ll ) do con- 
sider counts at these wavelengths from BLAST data. It is 
interesting to note that all five of these models over-predict 
the Herschel 250 — 350 /xm source counts (|piiver et al.l 
|2010| ; IClements et~all |2010|; iGlenn et al.l |2010| ). We find 



that the IValiante et alj (|2009f ) and IFranceschini et al] 
(|201dh models provide a decent fit (within 3a) to the 
1.1 mm source count s in all but the two lowest flux bins. 
IValiante et al.l (|2009l ) did not use the source counts at 
1.1mm to constrain their models, which at the time of 
their publication were l imited to data from r elatively small 
and/or shallow fields (|Laurent et al] 120051; Bertoldi et al 



2007; lPerera et alj|20 08). In comparison. IFranceschini et al 
( 20101 ) did use constraints on the 1.1mm counts from the 



AzTE C/JCMT survey of COSMOS j|Austermann et al] 
2009), and their model is in somewhat better agreement 
with our combined-field counts. However, compared to the 
other AzTEC fields shown in Fig. Q] the AzTEC/JCMT 
survey targeted a smaller and considerably over-dense 
region with in COSMOS, so it i s not too surprising that the 
model from IFranceschini et al . (2010) n ow over-pr edicts the 
counts at low flux densities. Marsden e t al.l (|201ll ) used the 
counts determined from the AzTEC/SHADES fields to fit 



their model, but it predates the correc tion to the AzTEC 
transfer function JDownes et al.l l201ll ); it now therefore 
under-predicts the combined-field counts at flux densities 
3mJy < Shoo < 5mJy, and over-predicts the counts in the 
lowest flux bin, by > 3a. However, this model is in better 
agreement with the count s at Si 100 < 3mJy than those of 
IValiante et~ail (|2009l ) and IFranceschini et all l|2010l ). 

We consider two phenomenological models that include 
recent co nstraints on the 250-5 00 /im counts from Her- 
schel. The lBethermin et al.l (|201ll ) model under-predicts the 
1.1 mm source counts by > 3a from 2mJy < Shoo < 6mJy. 
Those authors used the 1.1mm counts fr om the AzTEC 
surveys of GOODS-S ( |Scott et al] |2010D and SHADES 
l|Austermann et al]|2010l ) to constrain their model; however, 
the published counts from those surveys were systematically 
low due to the error in our transfer function e stimate. On 
the other hand. lRahmati fc van der Wera (l201ll ) did not use 
the 1.1mm counts from published AzTEC surveys to con- 
strain their model; however, they did compare their best-fit 
model to these data and showed that they are largely con- 
sistent, though their model over-predicts the 1.1mm counts 
at Shoo < 3 mjy. 

As demonstrated in Section [<01 we expect the extracted 
1.1 mm counts to be moderately biased due to confusion 
effects. Given the small statistical errors on the measured 
counts from these surveys, these biases should be taken into 
account when assessing the agreement between these data 
and various models. In Fig. [3] we include one-sided extended 
error bars (in grey) to account for the bias in our mea- 
surements due to confusion effects. These error bars encom- 
pass the 68.3 per cent confidence interval on the corrected 
counts, where we have used the results from the simulations 
in Section l"3.2l to estimate correction factors to the measured 
counts. While this is only an approximate correction for the 
effects of confusion, it is sufficient for our purpose of qual- 
itatively comparing our measurements to predictions from 
galaxy evolution models. As can be seen in Fig. [3] extend- 
ing the uncertainty on the counts naturally brings the data 
and models into somewhat better agreement; however, the 
discrepancies at Shoo > 4mJy remain large, and the main 
conclusions from the comparisons to various evolution mod- 
els above are unchanged. 



6 DISCUSSION 

Given that galaxy evolution models using semi-analytical 
methods provide some insight into the physical processes 
occurring within galaxies (albeit with several simplis- 
tic assumptions), it is interesting to compare these pre- 
dictions with our measured counts. The best-fit semi- 
analytical model to the combined 1.1mm counts is that of 
Gran ato et al.l (2004) . That model was able to reproduce the 
850 /im source counts and their redshift distribution mea- 
sured at that time, as well as the A'-band luminosity function 
of massive spheroids at z = 1.5. A key feature of their model 
is that feedback from SN is more effective in slowing down 
the rate of star formation in shallower potential wells, so that 
star formation progresses more rapidly within the most mas- 
sive halos. This scenario seems to be consistent with recent 
evide nce in favour of "downsizing" of SMGs fe.g lDve et al] 
I2OO81 ). where, contrary to the hierarchical collapse of dark 
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Figure 3. Comparison of the observed, Euclidean-normalised differential source counts from AzTEC surveys and predictions from galaxy 
evolution models. The differential counts derived from the combined six AzTEC blank fields are shown as circles (same as Fig. [TJ, where 
the black error bars show the 68.3 per cent confidence interval, including systematic uncertainties from cosmic variance (Section I3.3II 
and flux calibration (Section I3.4II . The one-sided extended error bars shown in grey encompass the 68.3 per cent confidence intervals 
including corrections to the measured counts due to bias from confusion effects, as discussed in Section [5] The squares show the bright 
source counts deriv ed from ACES fields, as described in Section [4] The triangles show the 1.4 mm source counts from the SPT survey 
IjVieira et al]|2010h . excluding nearby IRAS galaxies and sources whose SEDs are dominated by synchrotron emission. The SPT data 
have been scaled to 1.1mm assuming a spectral index of a = 2.65. The curves correspond to predictions from various semi-analytical 
and phenomenological models taken from the literature, as listed in the legend. 



matter, star formation in the early Universe predominately 
takes place within the most massive syst ems and progresses 
to lower mass systems at later times (e.g.lCowie et al.lll99£ 



Bundy et al.ll200rj; iFranceschini et al.ll200rj ; iMobasher et al 



.a 

06; 



20091 ; iMagliocchetti et a.l.ll201lT >. This idea is supported by 
the strong evolution of the luminosity function determined 
from phenomenological models, which implies that the most 
luminous star-forming galaxies (with IR luminosities Lib. > 
10 11 Lq) dominate at z > 1.5, while normal galaxies domi- 
nate at lower redshifts. 

It is interesting to note that, of the phenomenologica l 
models discussed here, only that of lFranceschini et al.l (|2010f ) 
is consistent with the turnover in the Euclidean-normalised 
counts at Shoo J; 2mJy; however, like all of these mod- 
els, it significantly over-predicts the counts at these faint 
flux densities, which have until now been only poorly con- 



strained. This observed turnover in our data is statistically 
significant. If we fit our measured counts near the appar- 
ent peak at Snoo = 2.4 and 3.4 mjy assuming no evolu- 
tion (i.e. flat in Euclidean- normalised space), we find that 
the observed counts at Siioo = 1.4 mjy fall short of this 
no-evolution model by 5a. The discrepancies between the 
predictions from these phenomenological models and the 
observed counts at faint flux densities may be highlight- 
ing our limited knowledge of both the SEDs and the den- 
sity of low luminosity (Lfir < 10 10 L©) galaxies in the lo- 
cal Universe and, i n turn, the faint end of the local lumi- 
nosity function f see IChapin. Hughes &: Aretx aga 2009, and 
references therein). In their phenome nological model l ing o f 
the 70-1100 fira detected populations. iMarsden et al.l (|2011h 
find that the model which best fits the observed counts and 
redshift distributions over-predicts the CIRB at these wave- 
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lengths. They explore whether simple modifications to the 
SEDs of low luminosity galaxies - to which the source counts 
are not sensitive, but which dominate the CIRB - can im- 
prove the fit to the CIRB. They find that by assigning the 
warmest SEDs to Lfir < 10 9 L0 galaxies, they can bring the 
CIRB prediction into agreement with the observed value. It 
is therefore possible that poor knowledge of the SEDs of 
local faint galaxies also limits how well phenomenological 
models can predict the faint end of the 1.1 mm source counts. 
Alternatively, evolutionary models that consider two dis- 
tinct populations evolving separately - where low luminosity 
galaxies evolve less strongly and thus reduce the num ber of 
cold galaxies at high redshift (e.g. lValiante et al.ll2009l ) - can 
also successfully match both the observed source counts and 
the CIRB. If this is the case, the turn over at the faint end 
of the AzTEC 1.1mm counts may be providing important 
new information on the evolution of low luminosity systems. 
At the other extreme, accounting for the observed num- 
ber density of bright SMGs poses a significant challeng e for 
existing theoretical models fe.g. lHavward et al.lbOllT ). As 
discussed in Section [H there is evidence for a bias in the 
observed source counts at Shoo > 15 mjy arising from back- 
grou nd SMGs that are strongly lensed by foreground struc- 
ture l|Vieira et al.ll2010f l. Although lensing is not expected to 
significantly bias the observed 1.1 mm counts at lower flux 
densities, current modelling of this effect is necessarily sim- 
plisti c and not well informed by obser v ations (Negrcll o et al.l 
20071 : IPaciea. Scott & . Chapir] 120081: [Lima. Jan fc Devlin! 



20101; iBetherminet al.ll201lP ). lAustermann et all J20ld) and 
Aretxaga et al.l ( 201ll ) have shown that Shoo > 5mJy 
sources are spatially correlated with z < 1.1 optical/IR 
galaxies in COSMOS, suggesting that galaxy-galaxy and 
galaxy-group lensing at moderate amplification levels may 
bias t he observed c ounts high, even at more modest flux den- 
sities. |Wanj2^et_alJ ([2011a) find a spatial correlation between 
foreground optical/IR galaxies from the Sloan Digitized Sky 
Survey (SDSS) and Spitzer/IRAC with high-redshift SMGs 
detected by Herschel in the Lockman-SWIRE field; since the 
redshift distributions of these different populations do not 
overlap significantly, this is strong evidence that the correla- 
tion arises from gravitat i onal le nsing. We have repeated the 
analysis of IWang et al.1 lJ2011al ) for the individual AzTEC 
fields for which SDSS and/or IRAC data are available, and 
find that only the COSMOS field shows a tentative (~ 2<r) 
correlation with the l ow-redshift galaxy samp les - consistent 
with the findings in lAretxaga et al.l (|201lh . Since we are 
sampling much smaller fields than the Herschel Lockman- 
SWIRE, our statistical power is limited. Still, these results 
hint that lensing may have a significant effect on the ob- 
served source counts, even at moderate flux densities. 

The observed counts can also be significantly biased 
if a large fraction of the SMGs detected as single point 
sources in low-resolution surveys, such as those carried out 
with AzTEC, LABOCA, and Herschel, are actually mul- 
tiple systems blended by the large beam. High-resolution 
(~ 2 arcsec) interferometric imaging with the Submillimeter 
Array (SMA) of an unbiased, flux-limited (Shoo > 5.5 mjy) 
sample of 15 AzTEC-detected SMGs discovered in COS- 
MOS showed that only two (13 per ce nt) are resolved into 
two, physically unassociated galaxies (| Younger et al.ll2007l . 
2009). On the other hand, simulations designed to match the 
observed source counts suggest that the fraction of multiple, 



blended galaxies for SMGs detected in low-reso lution obser- 
vations can be even higher. IWang et al.l 1)20110 ) predict that 
~ 1/3 of the Sg50 > 5 mjy SCUBA sources are actually mul- 
tiple galaxies blended by the beam. Similarly, usi ng the same 
kind of simulations described in Section 13.21 IScott et al.l 
(|20101 ) demonstrated that ~ 25 per cent of SMGs detected 
as single point sources in the confusion-limited AzTEC map 
of GOODS-S are likely to be two or more sources blended to- 
gether. However, our simulations in Section T3. 21 suggest that 
this has only a small effect on the measured source counts. 
There are nevertheless potential biases to the observed 
source counts that we have not studied in our simulations: 
in particular, the effects of galaxy clustering. There are 
currently only weak c onstraints on the c lustering strength 



of bright SMGs (e . g. iBlain et aT] 



iLindner et~aI1l201ll ; William s et al 



2004 IWeifi et al.1 120091 : 
20111 ). and the cluster- 



ing strength of faint mm-selected sources below the confu- 
sion limit has not been measured. If these faint, but much 
more numerous galaxies are strongly clustered, the biases 
to the measured 1.1 mm source counts could be much larger 
than those predicted by our simulations in Section f3. 2 1 Our 
lack of knowledge on the clustering properties of sub-mjy 
SMGs precludes a more rigorous study of the potential bi- 
ases to the measured source counts from galaxy clustering. 
The fraction of blended SMGs detected in low-resolution 
surveys, and whether this leads to significant biases to the 
source counts, can be addressed in future observations with 
ALMA. 



7 SUMMARY AND CONCLUSIONS 

We have combined previously published data from six blank- 
field surveys at 1.1 mm taken with AzTEC, totalling 1.6 deg 2 
in area with rms depths of 0.4-1.7 mjy beam -1 , in order 
to derive the strongest constraints to date on the 1.1mm 
extragalactic source counts from Shoo = 1 — 12 mjy. We use 
the well tested bootstrap sampling method on the source 
catalogues to derive the counts, which allows for an accurate 
estimate of statistical errors, including correlations among 
the selected flux bins. Given the large total area sampled, 
we expect a systematic uncertainty in the counts arising 
from cosmic variance of only 3.9 per cent. As discussed in 
Sections 13.31 and 13.41 it is important to include systematic 
uncertainties from cosmic variance and flux calibration (~ 
5 per cent) into the total error on the observed source counts, 
since these must be considered when using the counts to 
co nstrain parameters in g alaxy evolution modelling (as done 
in iBethermin et al.l l201ll for example) . We have included 
these systematic uncertainties in the total errors reported 
on our combined-field source counts by incorporating them 
directly into the bootstrap sampling method, and we list 
these in Tables H El andg] 

Comparing the observed 1.1mm source counts to pre- 
dictions from several galaxy evolution models, we find that 



the agreement at flux densities Sn 



4mJy is generally 



good. Given that most of these models had been fit to bright 
source counts at 850 fim and/or 1.1 mm from previously pub- 
lished surveys, this agreement is expected. However, we find 
significant (> 3<r) discrepancies between the combined-field 
1.1 mm counts and many of these models at Shoo i$ 4mJy. 
Similarly, with the exception of the most recent phenomeno- 
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logical models that include constraints from Herschel sur- 
veys, many of the models over-predict the counts at 250- 
500 j-im. The data presented here provide strong constraints 

- highly complementary to those from Herschel surveys at 
shorter wavelengths - that should be used in future mod- 
elling of the formation and evolution of IR-bright galaxies. 

Of the s emi-analytical models considered in this paper, 
the model of iGranato et al.l (J200J) provides the best match 
to the 1.1 mm source counts. This model is consistent with 
downsizing, in which the bulk of star formation activity pro- 
gresses from more massive to less massive galaxies over time 

- a scenario that is also supported by most phenomenolog- 
ical models that can describe the counts and redshift dis- 
tribution of I R-bri ght sources. However, for the model of 
iBaueh et al.l l|2005h . where the build-up of stellar systems 
is consistent with hierarchical formation, only minor mod- 
ifications are required to improve its agreement with the 
observed source counts at 1.1 mm. A better understanding 
of the physical processes of gas cooling and feedback are 
needed in order to provide information on which scenario 
bests describes the assembly of massive galaxies. 

For the first time, we have been able to strongly con- 
strain the 1.1 mm counts at 51100 < 3mJy, and we measure 
a turnover in the Euclidean-normalised counts at 5*iioo < 
2mJy which none of the evolutionary models considered 
here are able to reproduce. This either reflects our limited 
knowledge of the SEDs of low luminosity galaxies in the 
local Universe, or motivates modifications to the evolution 
of faint galaxies at high redshift. Wi de area surveys at 60- 
500 /im with Herschel (e.g. H-ATLAS. feales et al.l2010t ) and 
future surveys at 450 /im with SCUBA-2 and at 1.1 mm with 
AzTEC on the Large Millimeter Telescope (LMT) will pro- 
vide measurements of the SEDs of a large, unbiased sample 
of faint nearby galaxies, and in turn, allow for improved 
modelling of the evolution of these systems at high redshift. 

While there is already considerable evidence that the 
counts at 5*1100 > 15mJy are biased high by strong lensing 
effects, some groups have also demonstrated that galaxy- 
galaxy and galaxy-group lensing with moderate amplifica- 
tion may affect the observed counts at more modest flux 
densities (Shoo ~ 5mJy). Furthermore, galaxies detected 
as single point sources in low-resolution surveys, such as 
those taken with AzTEC, can also bias our measurements of 
the source counts if a significant fraction of these "sources" 
are really two or more galaxies blended by the large beam. 
The effects of lensing and blending on the measured source 
counts are poorly understood, and may lead to inaccurate 
predictions from galaxy evolution models. High-resolution 
imaging of statistically significant samples of SMGs with 
ALMA will be possible in the near future, and these will 
greatly aid in quantifying the degree to which lensing and 
blended galaxies bias the observed source counts. 

Compared to observations at shorter submm wave- 
lengths from Herschel, surveys at 1.1mm sample (on aver- 
age) galaxies at higher redshifts, as the negative fc-correction 
for observations at longer wavelengths extends to z ~ 10. In- 
deed, a growing number of 1.1 mm -bright SMGs are found 
to be at z > 4 (iDaddi et all I2009allb1 ; ICoppin et al.) [200p ; 



iRiechers et al.ll2010l ; ICapak et al.ll201ll ; ISmolcic et alfeoilT ) 
The source counts at 1.1 mm presented here, and those de- 
rived from future surveys with AzTEC on the LMT, will 
therefore provide crucial information on the evolution of 



star-forming galaxies within the first ~ 2 Gyr after the Big 
Bang, where observations at shorter submm wavelengths 
provide few constraints. 
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